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ABSTRACT 

We use high-resolution N-body simulations to develop a new, flexible, empirical ap- 
proach for measuring the growth rate from redshift-space distortions (RSD) in the 
2-point galaxy correlation function. We quantify the systematic error in measuring 
the growth rate in a 1 Gpc^ volume over a range of redshifts, from the dark 
matter particle distribution and a range of halo-mass catalogues with a number den- 
sity comparable to the latest large-volume galaxy surveys such as the WiggleZ Dark 
Energy Survey and the Baryon Oscillation Spectroscopic Survey (BOSS). Our sim- 
ulations allow us to span halo masses with bias factors ranging from unity (probed 
by emission-line galaxies) to more massive haloes hosting Luminous Red Galaxies. 
We show that the measured growth rate is sensitive to the model adopted for the 
small-scale real-space correlation function, and in particular that the "standard" as- 
sumption of a power-law correlation function can result in a significant systematic 
error in the growth rate determination. We introduce a new, empirical fitting func- 
tion that produces results with a lower (5-10%) amplitude of systematic error. We 
also introduce a new technique which permits the galaxy pairwise velocity distribu- 
tion, the quantity which drives the non-linear growth of structure, to be measured 
as a non-parametric stepwise function. Our (model-independent) results agree well 
with an exponential pairwise velocity distribution, expected from theoretical consid- 
erations, and are consistent with direct measurements of halo velocity differences from 
the parent catalogues. In a companion paper we present the application of our new 
methodology to the WiggleZ Survey dataset. 

Key words: large-scale structure of Universe, cosmological parameters, cosmology: 
theory 



1 INTRODUCTION 



The growth rate of cosmic structure is a key parameter 
which quantifies the cosmological model. While distance- 



redshift probes such as Type la supernovae (e.g. Riess et al 



19981 IPerlmutter et al.| |1999| 


Kowalski et al.| |2008| |Hicken| 


et al. 


2009 


AmanuUah et al. 


20101 or baryon acoustic os- 



2om 



cillations (e.g. Blake & Glazcbrook, 2003; Sco & Eiscnsteinl 
Eisenstein et al.[ |2005, ,Percival et al., ,2010, ,BlaE| 



et al. 2011 1 give us information about the cosmic expansion 
history, and observations of the Cosmic Microwave Back- 
ground radiation allow the study of the composition and 



physics of the early Universe (e.g. Komatsu et al. 20111, 
the growth rate describes how the small density perturba- 
tions present in the early Universe evolve to form the large- 
scale structure which populates the Universe today. Hence 



the growth rate provides a fundamental test of the laws of 
gravity which operate in the expanding cosmos. 

Measurements of the growth rate have assumed a spe- 
cial importance as evidence has accumulated that the ex- 
pansion of the Universe has entered a phase of acceleration. 
Two principle explanations have been put forward for accel- 
erating expansion. The first is the presence of some unknown 
material constituent of the Universe with the exotic property 
of negative pressure, whose energy density has become dom- 
inant within the last half of the age of the Universe. This ma- 
terial is known as "Dark Energy" . The second explanation is 
that our current theory of gravity. General Relativity, must 
be modified on large scales to account for the observations 
without invoking exotic constituents. Various methods for 
modifying General Relativity have been explored (see |Tsu-| 



jikawa 2010 1. Measurements of the growth rate of structure 
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over different cosmic epochs can help discriminate between 



these two physical interpretations of the observations (Lin- 
I der fc Jenkins, "2^3] [Linder fc Cahnl [20071 |Guzzo et al"] 
2008, Wang , 2008| ). 



The growth rate of structure within a given cosmologi- 
cal model can be derived by solving the differential equation 
for the linear density of matter perturbations S at scale fac- 
tor a in an expanding Universe ( Peebles 1980") . The growth 
rate at redshift z is defined by f{z) = dlnS/dlna and is 
well-approximated in a variety of dark energy models by 
f{z) = Q,m{z)'' , where ^lm{z) is the matter density relative 
to the critical density, and 7 is a phenomenological parame- 



ter which takes the value 0.55 for General Relativity (Linder 



& Cahn 20071 



A powerful method for measuring the growth rate is to 
exploit the anisotropic signature it imprints in the clustering 
within galaxy redshift surveys, known as redshift-space dis- 



tortious (RSD, see e.g. Hawkins et al. 2003; Tcgmark et al. 


2004| 2006| 


da Angela et al. 2008| Guzzo et al.,,2008 Oku-| 


mura et al. 


|2008| Cabre & Gaztanaga| 2009| Blake et al.| 


2011) . The 


growth of structure is driven by coherent flows 



of matter into clusters and superclusters, and the result- 
ing coherent galaxy peculiar velocities produce a correlated 
signature in the pattern of 2-point galaxy clustering. The 
galaxy clustering anisotropy on large scales is parameter- 
ized by a measured distortion parameter /3, which is related 
to the growth rate by / = 6/3, where b is the bias parameter 
which relates galaxy overdensity to matter overdensity. 

An alternative method to map out the clustering pat- 
tern and measure the growth rate is weak gravitational lens- 
ing, which involves modelling the observed correlated align- 
ment of distant background galaxy shapes by the foreground 
large-scale structure. However, the current sensitivity of the 



method is limited (e.g. Bernstein 20091, despite extremely 



promising future prospects. Other methods that have been 
developed to constrain the growth rate include the luminos- 
ity function and gas-mass fraction of X-ray selected clus- 



ters (Rapetti et al. 20091, and galaxy bulk flows measured 



in the local neighborh ood ([Abate fc Lahav[ [2008| [Watkins 



et al. 2009 Nusser & Davis 2011 1. However, at the present 



time redshift-space distortions in galaxy clustering represent 
the most accurate method for measuring the growth rate of 
structure. 

Utilizing redshift-space distortions as a cosmological 
probe has two requirements. Firstly, a large galaxy redshift 
survey must be performed in order to reduce the sources of 
error in the measurement, cosmic variance and galaxy shot 
noise. Secondly, reliable models must be fltted to the clus- 
tering measurements. The second issue, correlation-function 
modelling, is the subject of the current paper, whereas ap- 
plications of these techniques to the WiggleZ Dark Energy 



Survey ( Drinkwater et al. 2010 1 are presented in a compan- 



ion paper. 

In order to model the two-dimensional redshift-space 
galaxy correlation function transverse and parallel to the 
line-of-sight, fs((T, tt), we must model the pattern of pecu- 
liar velocities acquired by galaxies as part of the growth of 
structure. There are (crudely speaking) two sources for these 
peculiar velocities: the random velocity a galaxy possesses 
with respect to its own group or cluster, and the velocity the 
galaxy acquires as part of bulk flows into bigger structures 
such as superclusters. These two effects distort the shape 



of the measured (o", tt) in a characteristic way which can 
be modelled and split from the underlying isotropic real- 
space correlation function. The challenge is to measure the 
redshift-space distortion parameter due to bulk flows, /3, in 
a manner free from systematic error, and relate that to the 
growth rate of structure using parallel measurements of the 
galaxy bias factor. 

The aim of this paper is to address the obstacles which 
must be overcome to obtain reliable measurements of the 
cosmic growth rate from modelling the correlation function, 
namely: (1) provision of an accurate model for the under- 
lying isotropic real-space correlation function, (2) modelling 
the non-linear and random effects of galaxy pairwise veloci- 
ties, (3) determining an accurate measurement of the galaxy 
bias factor. We assume throughout that the background 
cosmological parameters are known, for example from Cos- 
mic Microwave Background observations, hence we neglect 
Alcock-Paczynski distortions in this study. 

There currently exists no model which is able to de- 
scribe accurately the matter distribution and its clustering 
properties over the complete range of scales of our interest: 
— 50h~^ Mpc. At least three different physical regimes can 
be identified which shape the galaxy correlation function 
and redshift-space distortions in this range: the non-linear, 
quasi- linear and linear regimes. On large scales (> 20 
Mpc), perturbation theory provides a good description of the 



growth of the clustering pattern | 


3ernardeau et al. 


2002 


Crocce fc Scoccimarroj |2006| 


Nishimichi et al.||2009 


Carl-j 


son et al. 


2009 


Taruya et al. 


200c 


) 1 , a continuity equation 



holds between velocity and density, and the redshift-space 
distortion pattern is well-described by a simple anisotropy 



known as the "Kaiser limit" (Kaiser 19871. At intermediate 
scales (< 20 Mpc), in the quasi-linear regime, pertur- 
bation theory breaks down (Peebles 1980 Hatton fc Cole[ 



1998 Ratcliffe et al. 1998 Landy 2002 Scoccimarro 2004 



and significant corrections are needed to the Kaiser limit 
formulation. At small scales (< 3h~^ Mpc), in the non- 
linear or "1-halo" regime, the matter distribution is virial- 
ized in groups and clusters of galaxies and high velocity dis- 
persions imprint the "fingers-of-god" feature into the corre- 
lation function. All of these physical ingredients leave differ- 
ent signatures in the correlation function, currently making 
it difficult to construct a physically-motivated model which 
is valid across the entire range of scales. In this study we 
instead focus on empirical approaches for modelling the ob- 
servations. 

N-body dark matter simulations, and the halo cata- 
logues that can be constructed from them, are a powerful 
method for modelling the full range of non-linear processes 
described above. Because the cosmological parameters and 
input growth rate in the simulations are known, we can 
study the systematic errors that arise from using particular 
algorithms to extract these observables. Simulations allow 
us to see in detail pieces of information that are typically 
hidden in the real data samples, such as the velocity distri- 
butions in virialized clusters or the shape of the real-space 
correlation function. 

Other studies in the literature have used simulations 
to test and improve models fitted to RSD. In |Cabre| 



fc Gaztaiiaga (2009) the pairwise velocity distribution in 
redshift-space was studied in detail and found to be scale- 
dependent, leading to the inclusion in the model of two inde- 
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pendent velocity dispersion terms a pplying at scales smaller 
and larger than 2 Mpc (see also 



Slosar et al. 



20061 



Jen- 



nings et al. ( 2011 1 used N-body simulations to test models of 
varying complexity for recovering the true value of RSD for 
different cosmologies via the power spectrum. They demon- 
strate that linear models by themselves do not extract an un- 
biased growth rate (see also Matsubara 2008 Taruya et al. 



20091 


Kwan et al.| 


et al. 


2012 


1. 



Kwan et alj [20TT] [Okumura fc Jing[ |2011[ |Bianchi| 



In our work we use a new set of N-body simulations the 
Gigaparsec WiggleZ (GiggleZ) Simulations (Poole et al., in 
prep) - to develop a new, empirical method for extracting 
the growth rate of structure from a range of halo catalogues 
and redshifts. The distinguishing feature of these simulations 
is their low particle mass in comparison with most large- 
volume simulations, which is appropriate for modelling the 
relatively low-mass haloes probed by emission-line galaxies 
mapped by the WiggleZ Dark Energy Survey. In contrast to 
the studies cited above, which largely considered the galaxy 
power spectrum, we base our analysis on the galaxy corre- 
lation function. This is a commonly-used statistic for quan- 
tifying galaxy clustering which has several merits including 
that (1) different physical processes (such as shot-noise) are 
confined to distinct sets of scales, and (2) it is less sensitive 
than the power spectrum to modelling the survey selection 
function. In the absence of a complete physical model, it has 
been standard in the literature to fit a power-law model as 
the real-space correlation function when modelling the data 
in redshift-space, and include non-linearities via an expo- 
nential pairwise velocity distribution function. In this paper 
we use the GiggleZ simulations to critically examine these 
assumptions, proposing a new, improved empirical fitting 
function and new techniques for studying the non-linear ve- 
locity dispersion. 

Our paper is structured as follows: in Section [2] we de- 
scribe the N-body simulations we employ in more detail. 
In Section [3] we outline how the 2-point correlation func- 
tion is measured from the simulation data, and in Section 
|4] we specify the models we fit to these measurements. In 
Sections |5] and |6] we describe the performance and results 
of these fits to the correlation functions measured from the 
dark matter distribution and halo-mass catalogues, respec- 
tively. In Section |6.2| we discuss the determination of the 
galaxy bias parameter from these data, and in Section [7] we 
summarize and discuss our findings. 



2 THE GIGGLEZ SIMULATION 

We conducted our analysis using the Gigaparsec WiggleZ 
Survey simulations (GiggleZ). The main simulation, which 
we utilize here, is a 2160^ particle dark matter simulation 
run in a periodic box lh~^ Gpc on a side. The result- 
ing particle mass of this simulation is 7.5 x 10^ h~^MQ 
which permits us to resolve bound systems with masses 
> 1.5 X 10^^ h~^MQ, facilitating studies of haloes with clus- 
tering bias factors ranging from near unity (e.g. galaxies in 
the WiggleZ Dark Energy Survey) to in excess of 2 (e.g. Lu- 
minous Red Galaxies (LRGs) in the Sloan Digital Sky Sur- 
vey). A WMAP-5 cosmology with {Qa,^ai ,^b, h,as,n) — 
(0.727,0.273,0.0456,0.705,0.812,0.960) was assumed for this 
simulation, with the initial conditions constructed to yield 
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Figure 1. The spread of maximum circular velocities and halo 
masses within the 6 GiggleZ simulation halo-mass catalogues an- 
alyzed in this study, for the 2 = snapshot. Each catalogue 
contains 250,000 position and velocity entries in a lh~^ Gpc 
box, generating a similar number density to that observed in the 
WiggleZ Dark Energy Survey or Baryon Oscillation Spectroscopic 
Survey. 



a CAMB ( [Lewis et al.[ |2000[ ) power spectrum for a start- 
ing redshift of = 49 using the Zeldovich approximation 



(Zel'dovich 1970 Buchert 19921 



Bound structures were identified using Subfind 
( [Springel et al. 2001 1, which uses a friends-of- friends (FoF) 
scheme followed by a substructure analysis to identify bound 
overdensities within each FoF halo. We use the Subfind sub- 
structures for all our analysis in this paper and use the value 
of each halo's maximum circular velocity Vmax as a proxy 
for mass. This choice of mass proxy was made to increase 
the robustness and reproducibility of our results, since it 
avoids many numerical effects and biases associated with 
specific halo finding schemes and halo mass definitions. We 
use the centre of mass velocities of each halo when comput- 
ing redshift-space distortions. 

In order to explore clustering systematics as a func- 
tion of halo mass, we rank-ordered the GiggleZ substruc- 
tures by their maximum circular velocities and selected con- 
tiguous groupings of 250,000 systems. From this we chose 
a series of 6 halo groupings, ranging in halo mass from 
WiggleZ galaxies to Sloan LRGs (with median values of 
Vmaa: = 130,160,190,230,260,300 km/s) with a number 
density 2.5 x 10"'* Mpc""^, which is well-matched to that 
of the WiggleZ survey and ongoing Baryon Oscillation Spec- 
troscopic Survey. Figure[T]illustrates the ranges of maximum 
circular velocities and halo masses contained in these six 
catalogues for the z — snapshot. We note that the com- 
pleteness limit of the halo catalogue is around 120 km/s at 
z = 0. 
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3 MEASUREMENTS OF THE 2D 

CORRELATION FUNCTION IN THE 
GIGGLEZ SIMULATION 

In this paper we quantify clustering using the two- 
dimensional 2-point correlation function. This statistic is de- 
termined by counting the number of unique galaxy pairs as 
a function of transverse and parallel separation to the line- 
of-sight (cr and vr), and comparing the result to a similar 
pair-count performed on a randomly-distributed catalogue. 
We averaged over 30 random catalogues, each containing an 
equal number of particles as the dataset. We used two esti- 
mators to measure the correlation function. The first is the 



Landy & Szalay ( 1993 1 estimator, the standard minimum- 



variance procedure performed with real galaxy catalogues 
and corresponding random catalogues. In the second estima- 
tion we exploit the fact that our simulation is a cube with 
periodic boundary conditions and only count data-pairs, es- 
timating the random pair-count using analytic methods. The 
results agree very closely, and the second method allows a 
much more rapid computation. We computed the pair counts 
in square (a, vr) bins with side 2 Mpc, up to 40 Mpc 
in a and 30 Mpc in tt (300 data bins). 

We obtained the data covariances using the jack-knife 
procedure, in which we divided the lh~^ Gpc cube into a 
number Njk of identical jack-knife regions. Singular Value 
Decomposition (SVD) analysis of the resulting covariance 
matrices is crucial for understanding their robustness in the 
fitting procedure. Figure [2] illustrates the spectrum of 
SVD eigenvalues for a range of choices of the number of jack- 
knife regions. We found that for our dataset, 7'^ jack-knife 
subvolumes produced noisy covariance matrices with both 
of the employed methods, and this effect is not completely 
ameliorated if the lowest eigenvalues are truncated in the re- 
constructed covariance matrix. SVD analysis showed that 
increasing the number of jack-knife sub-regions improves 
asymptotically the quality of the covariance matrix. In our 
case values in the range 10^ to 20^ are good choices. For 
our default choice Njk = 10^, we note that the size of each 
jack-knife region is 100 Mpc, significantly exceeding the 
clustering scales of interest. We checked that the growth- 
rate measurements presented in this paper do not depend 
significantly on whether we use the full covariance matrix, a 
truncated matrix in which the lowest-amplitude eigenvalues 
of an SVD decomposition are excluded, or a diagonal error 
matrix. 

We performed these measurements using simulation 
snapshots at redshifts z — and z = 0.6, in both real-space 
and redshift-space, for the dark matter distribution and for 
the 6 different halo-mass catalogues. We selected z = 0.6 
because this is the median redshift of the WiggleZ Dark En- 
ergy Survey, and z — 0.0 because here the non-linearity in 
the clustering pattern that we are modelling will be most 
significant. We randomly sub-sampled the dark matter cat- 
alogue for each snapshot to 10^ particles before performing 
the correlation function estimation. The expected values for 
the growth rate / for these 2 snapshots are 0.49 and 0.76 
respectively, corresponding to the ACDM cosmological pa- 
rameters used to construct the simulations. When generat- 
ing redshift-space positions we used a plane-parallel approx- 
imation, shifting co-ordinates along one axis. In Figure|3]we 
plot the 2D real-space and redshift-space correlation func- 
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Figure 2. The curve of eigenvalues in a Singular Value Decom- 
position of the covariance matrix which results as we vary the 
number of jack-knife subregions. The robustness of the covari- 
ance matrix improves asymptotically until we reach Njk = 20^, 
which is approximately the maximum number of sub-volumes in 
which the original data cube can be divided whilst still retaining 
dimensions bigger than the relevant scales in our 2-point correla- 
tion function measurement. The shape of these curves is similar 
for both methods of calculating the covariance matrix described 
in Section |3] In the case of lower numbers of jack-knife subdivi- 
sions, we were forced to truncate the lowest eigenvalues of the 
covariance matrix to attain stable fits, whilst for Njk > 10^ full 
and truncated covariances gave consistent results. 



tions of the GiggleZ simulation dark matter distribution for 
the z — 0.6 snapshot. In Figure|4]we plot the corresponding 
2D redshift-space correlation functions of the 6 halo-mass 
catalogues. 



4 MODELLING THE REDSHIFT-SPACE 
CORRELATION FUNCTION 

4.1 Constructing the model 

The fact that large-scale galaxy surveys observe redshifts, 
not distances, implies that their clustering pattern is dis- 
torted by galaxy peculiar velocities. These peculiar veloci- 
ties are generally modelled in the correlation function by a 
combination of the two effects which dominate in the large- 
scale and small-scale limit: the large-scale coherent fiow of 
galaxies into clusters and superclusters, and the small-scale 
random motions of galaxies within virialized structures. The 
large-scale effects of coherent flows on the power spectrum 
and correlation function can be described by the standard 



treatment of Kaiser (19871 and Hamilton (19921 and the 



small-scale random velocity distribution can be introduced 
by convolving with a function f{v), as summarized for ex- 
ample by Hawkins et al. (20031) 



We note some potential systematic errors in this ap- 
proach, which we discuss in turn in the remainder of this 
Section. Firstly, the isotropic real-space correlation func- 
tion 5, (s) must be modelled reliably in order to extract the 
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Figure 4. Measurements of the 2D redshift-space correlation function for 6 different GiggleZ simulation halo-mass catalogues from low 
mass (top left) to high mass (bottom right) in the z = 0.6 snapshot. 



anisotropic signature. Historically, a power-law has been em- 



ployed (Hawkins et al. 2003 Madgwick et al. 2003 Cabre 



& Gaztafiaga 20091. However, with increasing quality of 
data and simulations, a power-law has become a bad ap- 
proximation to the true non-linear clustering pattern. We 
discuss some improvements below. Secondly, the non-linear 
behaviour of small-scale peculiar velocities is entirely mod- 
elled by the single function ,f{v). However, in reality this 
function is describing a scale-dependent combination of at 
least two physical effects; the virialized motion of galaxies 
within haloes, and the non-linearity in the coherent flows of 
galaxies which damp the velocity power spectrum on quasi- 



linear scales ( |Sheth| [T996} [Slosar et al.[ |2006| ). In detail, a 
single function is unlikely to provide a good match in both 
regimes, and indeed more complex models have been con- 
sidered for matching small-scale data ( [Cabre &: Gaztafiaga[ 
|2009j). Another possibility is to reduce the impact of system- 
atic modelling errors by excluding data at small scales from 
the fit ( [Hawkins et al.[ |2003[ l , although this comes at the 
price of throwing away a fraction of the data and in conse- 
quence obtaining a statistically poorer measurement of the 
growth rate. 
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4.2 Models for the real-space correlation function 

4.2.1 Fitting formulae from CAMB and halofit 

For a given set of cosmological parameters, the matter power 
spectrum (hence correlation function) at recombination can 
be numerically calculated by solving Boltzmann's transport 
equation; a popular publicly-available code which provides 
this solution is CAMB ( [Lewis et ar||2000| . The effect of the 
non-linear growth of structure at an arbitrary redshift can 
be incorporated using the "halofit" recipe calibrated by N- 
body simulations ( [Smith et al.[[2003l ). In this model we can 
determine the galaxy correlation function by combining this 
with a linear galaxy bias parameter fo as a simple normaliza- 
tion. We refer to the non-linear real-space galaxy correlation 
function generated in this manner as the "CAMB model", 
and by combining this real-space correlation function with 
the redshift-space distortion parameters described above we 
can determine the corresponding redshift-space correlation 
function: 



^s{a,7v) = CAMBmodel [6, /3, /(t;)] 



(1) 



We explore below the dependence of the results on the fit- 
ted range of scales. We obtained the input CAMB power 
spectrum using the WMAP-7 best-fitting cosmological pa- 



rameters (Larson et al. 2011 1 which are consistent with the 



input parameters for the GiggleZ simulation. 



4.2.2 Power-law correlation function 

In previous studies the real-space galaxy correlation function 
at small scales has often been modelled with a power-law 



form = 



and the full model for the redshift-space 



correlation function can be written as 

^s(a, tt) = PowerLawModel [7, tq, /3, f{v)] (2) 

4.2.3 Quadratic correlation function ((,r,QCF) 

Improving data from both galaxy surveys and numerical 
dark matter simulations has demonstrated that the real- 
space galaxy correlation functio n deviates from a pow er-law 
at scales beyond ~ 15 Mpc (Hawkins et al. 20031. Thus 



in the case of high signal-to-noise data, a power-law model 
produces a poor fit to the correlation function. This moti- 
vates our definition of a new empirical fitting formula with 
greater flexibility than the CAMB and power-law models. 
We introduce here the quadratic correlation function ( QCF) 
model 



ir{r) 



-7 + qlogio(:^) 



(3) 



where q is the additional quadratic parameter. This is just 
a simple quadratic equation in logarithmic space: 



y = a^hx + qx 



X = logio('') 



Our full model in this scenario is then 
€,s{o, tt) = QCFModel [7, ro, go, f{v)] 



y = logio(?r) (4) 



(5) 



Although this model is not physically motivated, it produces 
an impressive fit over the wide range of scales 1 — 50 
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Figure 5. Fit of tlie power-law, CAMB and QCF models (dashed, 
dotted and continuous line-style respectively) to the real-space 
correlation function measured from the GiggleZ simulation dark- 
matter catalogue at redshift 2 = 0.6. The best-fitting values of 
the parameters and reduced chi-squared values Xr listed in the 
Figure. This plot shows that a Power Law model does not accu- 
rately describe the real-space correlation function, which causes 
systematic errors in the determination of the growth rate and 
galaxy bias parameters in this case. The CAMB and QCF mod- 
els provide a better match to the measured correlation function. 



Mpc to both a suite of non-linear matter correlation func- 
tions generated by CAMB, and (as we show below) to the 
real-space correlation function of a range of halo-mass cat- 
alogues from the GiggleZ simulation. The flexibility of the 
QCF model enables us to achieve flts to the redshift-space 
correlation function with lower systematic errors. We note 
that the QCF model, when applied to a halo correlation 
function, can also accommodate a scale-dependent bias fac- 
tor. 

Examples of these fits are given in Figure [5] and Ta- 
ble |2] In Figure |5| a comparison is shown between the ID 
real-space correlation function measured from the GiggleZ 
simulation at redshift z = 0.6, and fits of the three real- 
space correlation function models defined above. In Table 
[2] similar fits are shown for the 6 halo-mass catalogues. We 
note that the power-law model always produces a poor fit to 
the halo correlation function. The CAMB model describes 
the correlation function of low-mass haloes well, but breaks 
down at higher mass due to the effects of scale-dependent 
bias. The QCF model produces a good fit to the data in all 
cases. 



4.3 Models for the pairwise velocity distribution 

In previous studies, two choices have been considered for the 
pairwise velocity distribution f{v) describing the random 
small-scale motions of galaxies which result in the "fingers- 
of-god" observed in galaxy redshift surveys: a Gaussian or 



exponential distribution ( Ratcliffe et al. 1998 Landy & Sza- 
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lay 1998 Hawkins et al. 2003 1. We write these functions as 

1 



1 f -v^ 



2a 



exp 



(6) 



where a is the standard deviation of both distributions. Intu- 
itively, a Gaussian distribution can be thought of as resulting 
from virialized motions within a single dark matter halo, and 
an exponential distribution can be thought of as a sum of 



these motions across haloes of different mass I Sheth 1996 1 



The exponential distribution generally produces a better fit 



to data and simulations ( Hatton & Cole |1998 [Landy fc 



Szalayl |1998| Ratcliffe et al.[ |1998| |Landy[ |2002[ ). For ex- 
ample, [Landy ( |2002| ) demonstrated that the Fourier trans- 
form of the pairwise velocity distribution of three different 
galaxy surveys (LCRS, 2dFGRS and SDSS) is better fit by 
a Lorentzian profile, implying an exponential distribution in 
configuration space. 

In this study we consider two variations. Firstly we ex- 
plore both possibilities by fitting to our data the weighted 
combination 



f^{v)=xf,{v) + (l~x)fg{v) 



(7) 



thus adding a final parameter, x, to our correlation-function 
fits, which we required to lie in the range Q < x < 1. 

Secondly we considered expressing the pairwise velocity 
distribution as a general stepwise function 



/(«) 



for 



-1 <= \v\ < Vi 



1,2, ...iV 



(8) 



for a number of intervals A*' from no = up to some max- 
imum pairwise velocity vm = v-max- The distribution is 
normalized such that f{v)dv — 1, which implies that 

2 ttiivi — Vi-i) = 1. The quantities A'^ and Vmax are set 

by hand (depending on the quality of the data) by inspect- 
ing the solutions and requiring that f{v) should be generally 
positive and smoothly decreasing from v — to Vmax- For 
our simulation dataset we typically obtain robust solutions 
using TV in the range 6 — 9 and Vmax between 1500 — 2500 
km s^^. 

For each combination of /3 and the parameters describ- 
ing the real-space correlation function, we obtained the set 
of coefficients which minimized the chi-squared statistic 
between model and data by solving an N x N linear system 
of equations. As explored in more detail below, we obtain 
stable, physically sensible solutions provided that we include 
correlation function measurements at small transverse sepa- 
rations cr in our fitted range (which have most sensitivity to 
virialized motions). Using this technique we can determine 
the real underlying shape of the pairwise velocity distribu- 
tion, and test whether or not the exponential or Gaussian 
models indeed provide a good description. We give algebraic 
details of this calculation in Appendix A. 



5 MODEL FITS TO THE DARK MATTER 
CORRELATION FUNCTION 

5.1 Growth rate of the dark matter 

Firstly we explored how accurately the three real-space cor- 
relation function models we have defined (CAMB, QCF, 
power-law) described the real-space dark matter correlation 
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Bias 


[h-i Mpc] 


xVdof 


xVdof 


xVdof 


1.006 





1.246 


6.909 


9.470 


1.006 
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1.234 


1.687 


9.425 


1.006 
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1.192 


1.161 


8.491 


1.004 
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1.019 


0.931 


6.438 


1.000 


8 


0.674 


0.704 


4.102 


0.999 


10 


0.633 


0.669 


2.559 



Table 1. This Table explores the effect of excluding the small- 
scale data bins from the fitting of the real-space dark matter 
correlation function. r„iin is the minimum value of the total sep- 
aration (r = \/a^ + TT^) for the data bins included in the fit. We 
note that the CAMB model provides a good fit to the real-space 
dark matter correlation function, even in the small-scale regime. 
The discrepancy in normalization between the fitted CAMB cor- 
relation function and the simulation is less than 1%. 



function of the simulation, before the addition of redshift- 
space distortions. One example of these results (for z — 0) is 
displayed in Figure [5] The minimum values of chi-squared, 
which are listed in Table [l] demonstrate that the CAMB 
model provides the best fit and a stable value of the galaxy 
bias ~ 1 for a variety of fitting ranges, whilst the power-law 
model provides a poor fit to the data. The QCF model also 
yields a good fit for all scales > 4h~^ Mpc. 

Next we included the effects of redshift-space distortions 
in the dark matter correlation functions of the z = 0.0 and 
z = 0.6 snapshots, and fit the clustering models described in 
Section [4] Figure [6] displays how the measured parameters 
for the two snapshots depend on the minimum transverse 
separation fitted, (Tmin, for the three different real-space cor- 
relation function models we are considering. We fix the max- 
imum transverse separation fitted at (Jmax = 40 Mpc. 
We performed the fits using a Monte Carlo Markov Chain 
(MCMC) procedure, exploring the multi-dimensional space 
of variables of the models and obtaining their joint and in- 
dividual probability distributions. The horizontal lines indi- 
cate the input values of the simulation, determined from its 
fiducial cosmological parameters. We find that the CAMB 
and QCF models are both able to recover the input growth 
rate with low systematic error and reduced ~ li 
Omin > 2h~^ Mpc. Ou the other hand, the assumption of 
a power-law model produces a significant systematic error 
and bad fit. 



5.2 Pairwise velocity distribution of the dark 
matter 

In Figure[6]we show the fitted values of the pairwise velocity 
dispersion parameters x and a, which are fit for each model 
jointly with the other parameters. The results strongly favor 
an exponential rather than a Gaussian velocity distribution 
{x ~ 1). The systematic discrepancy between the values of 
a fitted in the QCF and CAMB models indicate that there 
is some cross-talk between a and the shape of the real-space 
correlation function on small scales. 

We now compare these measurements with a direct de- 
termination of the shape of the pairwise velocity distribution 
f{v), which is possible using our stepwise fitting method de- 
scribed in Section 14.31 This constitutes a further check for 
systematic errors in the models which describe this distri- 
bution. Wc fit the stepwise distribution in 7 velocity bins up 
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Figure 6. Fits for the growth rate, bias, the pairwisc velocity dispersion, the parameter x and the resulting reduced x^j from the GiggleZ 
redshift-space 2D dark matter correlation function, for different values of the minimum transverse scale included in the fit, amin- Top: 
results from the z = 0.6 snapshot of the simulation. Bottom: results from the z = 0.0 snapshot of the simulation. Measurements are 
shown for three different real-space correlation function models. The CAMB and QCF models (red circles and green triangles) show good 
agreement with the expected theoretical simulation growth rate (represented by the horizontal line) and rest of parameters, while the 
power-law model fit (blue squares) is strongly affected by systematics.The parameter x controls whether the pairwise velocity distribution 
is modelled as an exponential or a Gaussian, with the data clearly favoring an exponential form {x 1). 



to a maximum velocity of 2000 km s^^ (our results are not 
sensitive to these choices). We assumed a CAMB real-space 
matter correlation function in these fits. 



6 MODEL FITS TO THE HALO 
CORRELATION FUNCTIONS 

We now turn our attention to fitting models to the halo 
correlation functions, which represent galaxies in the simu- 
lation. We explore how accurately our empirical models can 
describe the halo clustering in redshift-space. 



The results are displayed in Figure [7| comparing the 
best-fitting stepwise distributions to both the best-fitting 
exponential and Gaussian models, and to a direct measure- 
ment of the pairwise velocity distribution from the cata- 
logues. In the top row, the stepwise velocity distribution is 
fitted to the correlation function of the dark matter sub- 
sample, for 3 different values of the minimum value of a 
of the data bins included in the fit. In the bottom row the 
pairwise velocity distribution is directly measured from the 
dark matter catalogues, this time varying the maximum to- 
tal separation {R = V''"^ + tt^) of the pairs considered for 
this measurement. We set a maximum value of R for this 
calculation because we only want to estimate the dispersion 
within a single bulk flow, not between different bulk flows. 
The distribution obtained from the stepwise fltting process 
is in excellent agreement with its direct measurement from 
the catalogue, demonstrating the capacity of the stepwise 
approach to recover the pairwise velocity distribution from 
the measured correlation function. The distributions we ob- 
tain are consistent with an exponential and not a Gaussian 
shape for /(«), and we recover the standard deviation of 
the pairwise velocity distribution a ~ 300 km s~^, inde- 
pendently of specifying a model for f{v). This agreement 
between the model-independent stepwise measurement and 
the assumed functional form increases our confldence in the 
reliability of the exponential small-scale redshift-space dis- 
tortion model. 



6.1 Redshift-space distortion parameters of the 
halo catalogues 

Table [2] compares models and measurements in real space, 
for the 6 halo-mass catalogues from the GiggleZ simula- 
tion considered in this study. The power-law model consis- 
tently fails to match the measured correlation function. The 
CAMB model produces a good description of the clustering 
of low-mass haloes, but breaks down for higher-mass cat- 
alogues due to scale-dependent halo bias. The QCF model 
is flexible enough to produce a good flt to the real-space 
correlation function of all the halo catalogues. 

We then fitted the redshift-space distortion models to 
the 2D redshift-space correlation function measurements for 
these 6 halo-mass catalogues. In Figure |8] we show the re- 
duced statistic corresponding to fitting our three RSD 
models to the range amin < a < 40 Mpc, < tt < 
30 h-^ Mpc, as a function of the choice of the minimum 
value amin- Our motivation for this analysis is that system- 
atic errors in the recovery of the redshift-space distortion 
parameter are likely to be most serious at the lowest values 
of a, where non-linear velocity effects such as the "fingers- 
of-god" are most significant (note however |Bianchi et al. 
(2012| ), who show in a similar analysis that increasing amin 
does not make systematics completely vanish). We find that 
the power-law real-space model is a bad fit to the data, 
whereas the CAMB and QCF models show mutual agree- 
ment and better values. The QCF model provides the 
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Figure 7. Top row: The best-fitting stepwise functions for tlie pairwise velocity distribution of tlie 2D dark matter redsliift-space 
correlation function, assuming a CAMB real-space model correlation function. The fits are repeated for 3 different values of the minimum 
transverse separation (Tmin- Bottom row: direct measurements of the pairwise velocity distribution in the dark matter catalogue, for 
different values of the maximum pair separation R. In both panels we fit Gaussian (green) and exponential (blue) functions to the 
histograms. In the upper panel we restrict the fitting to the positive and monotonic part of the data. From these fits we get the value 
of the standard deviations. This plot demonstrates that we are able to recover the pairwise velocity distribution correctly using this 
method, and we find excellent agreement between the measured velocity distribution function and the exponential model. 
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Table 2. This Table compares models and measurements in real 
space, for the 6 halo-mass catalogues from the GiggleZ simulation 
considered in this study. The power-law model consistently fails to 
match the measured correlation function. The CAMB model pro- 
duces a good description of the clustering of low-mass haloes, but 
breaks down for higher-mass catalogues due to scale-dependent 
halo bias. The QCF model is flexible enough to produce a good 
fit to the real-space correlation function of all the halo catalogues. 



best fit to the clustering pattern, particularly for high-mass 
haloes. 

The best-fit values for the pairwise velocity dispersion 
a depend on the real-space correlation function model in 
a similar way to the fits to the dark matter catalogues. For 



these halo catalogues the x parameter, which was introduced 
to determine the shape of the pairwise velocity distribution, 
does not favor the exponential function over the Gaussian 
function as strongly as in the case of dark matter. We spec- 
ulate that a halo catalogue with a narrow range of masses 
will also possess a more uniform dispersion of pairwise veloc- 
ities, producing a closer match to a Gaussian function than 
before. In order to convert these results into growth rate 
measurements we require an estimate of the galaxy bias, 
which we consider in the next Section. 



6.2 Bias factor of the halo catalogues 

Modelling the RSD in the 2-point correlation function gives 
us one of the two quantities, /3, that are necessary to deduce 
the growth rate f = fib. We now discuss and test possible 
methods to measure the bias factor b, which describes the 
clustering of galaxies (or dark matter haloes) relative to the 
underlying matter distribution. For example more massive 
haloes, which are sampled preferentially from more clustered 
regions of the Universe, will possess a higher bias factor b 
and a correspondingly lower value of /3 (resulting in a less 
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Figure 8. The value of the reduced statistic of the best-fitting redshift-space distortion model as a function of the minimum value 
of (T included in the fit. Results are shown for each of the 6 halo-mass catalogues at redshift z = 0.6, for the three different RSD models. 



flattened large-scale redshift-space 2D correlation function). 
These more massive haloes are sites of early galaxy forma- 
tion which today can be observed as luminous red galax- 



ies (Eisenstein et al. 2001 Cabre & Gaztanaga 20091; less 



meissive haloes preferentially host blue, star-forming galaxies 



IMadgwick et al. 2003 Blake et al. 20101 



The bias factor cannot be deduced directly from the 
measured 2-point correlation function without some other 
assumption (such as the underlying amplitude of the real- 
space matter correlation function), although we note that 
such a determination may be possible using the 3-point cor- 
relation func tion ([Verde et al.| |2002[ IGaztafiaga fc Scocci-| 



2005[|Marin[|2011[ ). However, the use of N-body sim- 
ulation catalogues allows us to compare the bias measure- 
ments resulting from fits of a real-space correlation function 
model with those directly determined from the data via 

b^ir)=^G{r)/^DM{r) (9) 

where ^G{r) and £,DM{r) are the real-space correlation func- 



tions of the halo (galaxy) catalogue and dark matter distri- 
bution, respectively. We assume hereafter bias is a constant 
in our fits, which is a good approximation on large scales. 

Table |3] shows the consistency of the bias factors re- 
sulting from fitting a CAMB model to the halo-mass cata- 
logues, with the direct measurements of the bias obtained 
by dividing the different halo-mass correlation functions by 
the dark matter correlation function and using equation |9] 
We divided the two correlation functions in each different 
separation bin, and averaged over each measurement as an 
independent estimate of the bias. We note that the two in- 
dependent measurements are broadly consistent, with the 
main discrepancy appearing for small-scales due to system- 
atic non-linear effects. This discrepancy can be ameliorated 
by restricting the fitting range to exclude small scales, which 
are polluted by the highest-amplitude systematic errors. 

The broad agreement between the bias factors fitted 
by the CAMB model to the 2D redshift-space correlation 
function, and those measured directly from the ID real-space 
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Table 3. Comparison between the fitted and measured bias fac- 
tors for the 6 halo-mass groups, for snapshots at redshifts z = 0.0 
and z = 0.6. In general we recover consistent values for the bias, 
but we note a marginally-significant systematic discrepancy for 
higher-mass haloes. 



correlation function of the simulation, gives us confidence to 
use the CAMB bias factors to measure the growth rate of 
our halo catalogues. 



6.3 Growth rate of the halo catalogues 

The resulting growth rate measurements for the z — 0.6 sim- 
ulation snapshot, combining the separate determinations of 
the redshift-space distortion parameter /3 and the galaxy 
bias, are plotted in Figure [9] This figure illustrates the am- 
plitude of the systematic error in measuring the growth rate 
as a function of halo mass, minimum transverse scale fit- 
ted, and model adopted for the real-space correlation func- 
tion. All models contain a systematic error A/ « 0.05 for 
fmin < 6h~^ Mpc. However, our fitting procedure, and 
the adoption of the QCF model, allows us to recover ac- 
curate measurements when excluding the small-scale region 
[cmin > 6h~^ Mpc). Similar results are obtained for the 
z = 0.0 snapshot. We note in general that the model fits to 
halo catalogues described in this Section contain somewhat 
higher systematic errors than the fits to the underlying dark 
matter distribution. 



7 CONCLUSIONS 

We have measured the 2D correlation function of the dark 
matter and halo-mass catalogues for 2 snapshots of the Gig- 
gleZ simulation at redshifts z = 0.0 and z = 0.6, and fitted 
different clustering models to recover the growth rate. We 
carefully study the effect on the results of the range of scales 
fitted and the halo mass, spanning halo bias factors between 
unity (hosting emission-line galaxies) and high-mass haloes 
traced by Luminous Red Galaxies. We list our conclusions 
as follows: 

• The commonly-used power-law model for the real-space 
correlation function produces a poor fit to the clustering pat- 
tern and a systematic error in the resulting growth rate. A 
real-space correlation function based on a non-linear CAMB 
model does better, but breaks down for high-mass haloes. 
We introduce a new empirical model, the quadratic corre- 
lation function (QCF model, which has one more degree of 
freedom than a power-law) , which provides a better descrip- 
tion of the real-space correlation function (particularly for 
high-mass haloes which possess significant scale-dependent 



bias) and produces a measurement of the input growth rate 
with a lower (5 — 10%) amplitude of systematic error. 

• We introduce a new technique which permits the mea- 
surement of the pairwise velocity distribution as a step- 
wise function from the redshift-space correlation function. 
The pairwise velocity distribution measured directly from 
the simulation catalogues is consistent with this model, and 
matches closely to the exponential function expected from 
theoretical considerations. 

• We have quantified the amplitude of systematic error 
in the measured growth rate from our 1 Gpc^ simula- 
tion as a function of halo mass, the model employed, and 
the minimum transverse scale amin fitted. We find that for 
o"min < 6h~^ Mpc, the systematic measurement error us- 
ing our procedure is A/ ~ 0.05. The adoption of the QCF 
model allows us to recover accurate measurements of the 
growth rate from halo catalogues when excluding the small 
scale region {amin > Qh~^ Mpc). This is consistent with the 
recent analysis of Bianchi et al.l ( 2012 1. 



• We note that the halo correlation function contains a 
higher level of systematic modelling errors than the dark 
matter correlation function, due to scale-dependent galaxy 
bias. Our modelling allows us to recover the growth rate 
from the dark matter particle catalogue with no detectable 
systematic error. 

We conclude that N-body simulations are an essential 
tool for testing and developing methods to measure the cos- 
mic growth rate, and that our empirical techniques should 
be useful for modelling correlation functions measured in 
the latest large- volume galaxy redshift surveys. A compan- 
ion paper will represent the application of these methods to 
the WiggleZ Dark Energy Survey dataset. 
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APPENDIX A: STEPWISE VELOCITY 
DISTRIBUTION METHOD 

As an alternative to assuming a Gaussian or exponential 
pairwise velocity distribution, we explored fitting a stepwise 
function (/(f) = in a range of velocities) to the data, 
together with the other model parameters. In order to avoid 
introducing a large number of additional free parameters, 
we have developed a method which allows these fits to be 
performed in a relatively simple and fast way. 

The standard expressions for the correlation function 



multipoles in the large-scale limit (Kaiser 1987 Hamilton 
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Figure 9. The best-fitting value of tlie growtfi rate as a function of tlie minimum value of cr included in the fit. Results are shown for 
each of the 6 halo-mass catalogues for the z = 0.6 snapshot, for the three different RSD models. This figure shows we can recover the 
growth rate with minimal systematics when using QCF and CAMB models for the real-space correlation function. 
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Starting from the real-space correlation function, such as 
obtained from a CAMB matter power spectrum, we re- write 



equation |A1| in the form: 

Co(r) = Co{/3)Ma,7v) 
Ur) = C2(/3)J2(<T,7r) 
U{r) = C4(/3)J4(<T,7r) 



(A6) 



where the form of the functions d and Ji can be deduced 
by analogy with equation |A1[ and we solve the integrals 
numerically. Now replacing the velocity distribution f{v) by 
a stepwise function Ofc we obtain: 



/oo 
-oo 



AT 



^'((7, Tx')dv - 



5' (a, ■K')dv 



(A7) 



where we separate the negative and positive parts of the 
integral over v because of the loss of symmetry implied by 
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the relation 

v{z + l) 



TT = TT + 



Hiz) 



(A8) 



for the convolution with the pairwise velocity distribution. 
We now numerically calculate the following terms: 

Qo{cr,'K,k)= / Jo{a,n)Po{i^)dv 
Q2{cr,T:,k)= / J2{a,iv)P2{fJ,)dv 

Qi{a,'K,k)= I J4{<J,n)P4{id)dv (A9) 

where Pi{j-L) are the Legendre polynomials in terms of /i, 
the cosine of the angle between the bin position and the 
line-of-sight axis tt. We also define analogous functions Q' 
integrating over the negative part of the v axis. Now we 
combine these functions for a given value of P to create a 
series of stepwise models ^{a,n) = '}2'^kF{P)a,-n,k, where 
the final F„^,r,k terms correspond to: 

F(/3)a,,r,fe = Co(/3) [Qo(a,7r,fc) + Q(,(a,7r,fe)] 
+C2(/3) [Q2ia,-K,k) + Q2{a,-K,k)] 
+Ci{l3) [Qiia, TT, k) + Q'^ia, tt, k)] (AlO) 

Now we have a model which is linearly dependent on a^, 
which constitutes an N x N linear system of equations. Our 
equation is: 



-D. 



(All) 



where D is the data array, C is the covariance matrix, i and 
j represent bins in the correlation function data, and k is 
the index for the stepwise velocity distribution. 

The n*'' equation of the N x linear system is: 

^afc^Cr/ +F(/?),-fcF(;9)i,„] = 

Cr] [AF(/3),-,„ + AF(^)i,n] (A12) 

The linear system may be solved by conventional methods 
and the parameter space in /3 can be quickly explored in 

search of the minimum ^ . The normalization of the step- 
wise function gives in this case the factor between the fitted 
galaxy correlation function and the CAMB matter correla- 
tion function, i.e. the bias 6^. If the real-space correlation 
function is a power-law, the first set of numerical integra- 
tions may be replaced by analytical expressions, with the 
clustering length ro of the power-law absorbed into the nor- 
malization factor. 
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